version 18.0               // version control
set processors 8           // to ensure replicability across different numbers of cores
clear all                  // clear existing data
macro drop _all            // and macros, clean slate
set seed 20220909          // set seed

local pgm  "dst-Appendix_Table_A3b_Panel_C_sleep_deprivation_rd_epa_kernel" // file name
local who  "Muzhe Yang"                                                     // author
local dte  "2025-01-22"                                                     // created date
local dte2 "`c(current_date)'"                                              // last run date
local tag  "`pgm'.do, created by `who' on `dte', last run on `dte2'"

capture log close
log using "code\analysis\tables\appendix\\`pgm'.txt", text replace 
display "`tag'"

**# data prep ------------------------------------------------------------------------------------

use "data_clean\dst-data04_for_estimation_within_250_miles", clear
summ dist_to_border
tab wave, missing
preserve 
keep wave TractFIPS sleep_depr
reshape wide sleep_depr, i(TractFIPS) j(wave)
compare sleep_depr1 sleep_depr2  // two waves of data are identical, all from 2018
restore
keep if wave == 2
codebook TractFIPS

*## (1) create the 9 cells, following page 215 of the following paper:
*## Giuntella, O. and F. Mazzonna (2019). "Sunset Time and the Economic Effects of Social Jetlag: Evidence from US Time Zone Borders." Journal of Health Economics 65: 210-226.

assert !missing(centroid_lat)
assert !missing(region)
gen     cell = 1 if centroid_lat > 40                        & region == "eastern and central"
replace cell = 2 if (34 < centroid_lat & centroid_lat <= 40) & region == "eastern and central"
replace cell = 3 if centroid_lat <= 34                       & region == "eastern and central"
replace cell = 4 if centroid_lat > 40                        & region == "central and mountain"
replace cell = 5 if (34 < centroid_lat & centroid_lat <= 40) & region == "central and mountain"
replace cell = 6 if centroid_lat <= 34                       & region == "central and mountain"
replace cell = 7 if centroid_lat > 40                        & region == "mountain and pacific"
replace cell = 8 if (34 < centroid_lat & centroid_lat <= 40) & region == "mountain and pacific"
replace cell = 9 if centroid_lat <= 34                       & region == "mountain and pacific"
quietly tab cell, gen(cell_)
summ centroid_lat if cell == 1 | cell == 4 | cell == 7
summ centroid_lat if cell == 2 | cell == 5 | cell == 8
summ centroid_lat if cell == 3 | cell == 6 | cell == 9

*## (2) control variables

global x "pop white_pct black_pct hispanic_pct pop_under_18yr_pct pop_65yr_over_pct age_median educ_hs_pct educ_coll_pct married_pct hh_size hh_income_median home_value_median no_health_ins_pct unemploy_pct"

des dist_to_border dist_to_border_signed $x centroid_lat daylight_dur_max

**# estimation prep -----------------------------------------------------------------------

encode CountyFIPS, gen(county_fips)
egen nomiss = rowmiss($x dist_to_border centroid_lat daylight_dur_max)
assert !missing(StateAbbr)
local radius = 100
gen sample_selection = (nomiss == 0 & StateAbbr != "AZ" & dist_to_border <= `radius')

**# table ---------------------------------------------------------------------------------

rdrobust sleep_depr dist_to_border_signed if sample_selection == 1 & dist_to_border <= 50, ///
         c(0) covs($x centroid_lat daylight_dur_max cell_2-cell_9) kernel(epa) vce(cluster county_fips) all
estimates store m1
summ dist_to_border if e(sample)
scalar bw = e(h_l)
scalar kernel = e(kernel)
scalar donut_hole = "No"
quietly etable, replace ///
        column(index) ///
        cstat(_r_b,  nformat(%9.3f)) ///
        cstat(_r_se, nformat(%9.3f)) ///
        mstat(N,                   nformat(%9.0fc) label("Number of observations")) ///
        mstat(case_bw = bw,        nformat(%9.3f)  label("MSE-optimal bandwidth")) ///
        mstat(case_kernel = kernel,                label("Kernel function used")) ///
        mstat(case_donut_hole = donut_hole,        label("Use a 5-mile radius hole")) ///
        stars(0.10 "*" 0.05 "**" 0.01 "***", attach(_r_b) decreasing pvname("p-value") nformat(%9.2f)) showstars showstarsnote ///
        novarlabel nocenter

rdrobust sleep_depr dist_to_border_signed if sample_selection == 1 & (dist_to_border > 5 & !missing(dist_to_border)), ///
         c(0) covs($x centroid_lat daylight_dur_max cell_2-cell_9) kernel(epa) vce(cluster county_fips) all
estimates store m2
summ dist_to_border if e(sample)
scalar bw = e(h_l)
scalar kernel = e(kernel)
scalar donut_hole = "Yes"
quietly etable, append

rdrobust sleep_depr dist_to_border_signed if sample_selection == 1 & dist_to_border <= 50 & region == "eastern and central", ///
         c(0) covs($x centroid_lat daylight_dur_max cell_2-cell_9) kernel(epa) vce(cluster county_fips) all
estimates store m3
summ dist_to_border if e(sample)
scalar bw = e(h_l)
scalar kernel = e(kernel)
scalar donut_hole = "No"
quietly etable, append

rdrobust sleep_depr dist_to_border_signed if sample_selection == 1 & (dist_to_border > 5 & !missing(dist_to_border)) & region == "eastern and central", ///
         c(0) covs($x centroid_lat daylight_dur_max cell_2-cell_9) kernel(epa) vce(cluster county_fips) all
estimates store m4
summ dist_to_border if e(sample)
scalar bw = e(h_l)
scalar kernel = e(kernel)
scalar donut_hole = "Yes"
quietly etable, append

collect title "Table: The Impact of Circadian Misalignment Estimated by RD, Using the Epanechnikov Kernel Function" 

collect layout 
collect addtags top_row[1]
collect recode etable_estimates 1 = column1 2 = column1 ///
                                3 = column2 4 = column2 
collect layout (colname[Robust]#result[_r_b _r_se] result[N case_bw case_kernel case_donut_hole]) ///
               (top_row[1]#etable_estimates[column1 column2]#cmdset#stars)

collect label levels top_row ///
        1 "Treat = 1 for census tracts located east of the time zone border; Treat = 0 for census tracts located west of the time zone border", modify
collect label levels etable_estimates ///
        column1 "All time zones in the Contiguous United States" ///
        column2 "Eastern Time Zone and Central Time Zone", modify 
collect label levels cmdset ///
        1 "(1)" 2 "(2)" 3 "(3)" 4 "(4)", modify
collect label levels colname ///
        Robust "Bias-corrected RD estimates with robust variance estimator", modify

collect style header top_row, level(label)
collect style header etable_estimates, level(label)
collect style header cmdset, level(label)
collect style header colname, level(label)
collect style column, dups(center)
collect style cell result[_r_b N case_bw], sformat(%s)
collect style cell border_block[corner], border(top, pattern(double))
collect style cell border_block[column-header], border(top, pattern(double))
collect style cell border_block[row-header], border(bottom, pattern(double))
collect style cell border_block[item], border(bottom, pattern(double))

collect export "code\analysis\tables\appendix\\`pgm'.xlsx", replace

log close
exit